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We study the one-band Hubbard model at half filling with hopping asymmetry and its efi'ective 
model at finite but large U up to the second order of t/U. Two variational wave functions, the 
resonating valence bond (RVB) wave function and anti-ferromagnetic- RVB(AF-RVB) coexisted wave 
function, are studied by variational Monte Carlo method on LxL square lattices up to L=12. Based 
on these two wave functions, the phase diagrams for both models are presented. For the Hubbard 
model, we find that there is a metal-insulator transition when the hopping parameter tmix which 
changes the local double occupant vanishes while only a metal-insulator crossover is explored for any 
finite tmix- For the effective model in which the perturbation expansion is up to the second order 
of tmix/U , a clear metal-insulator transition can be identified for both variational wave functions 
and the phase diagram can be drawn accordingly. In both models, we find that the systems are 
dominated by AF-RVB wave function when U is large while RVB wave function is favored when U 
is small. 
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I. INTRODUCTION 

The one-band Hubbard model is a generic model for 
interacting electrons in the narrow-band and strongly 
correlated system^i. Especially, since the high temper- 
ature superconductivity was discovered in the cuprates, 
the Hubbard model on two-dimensional lattices as well 
as its strong coupling limit modeP, the t-J model, have 
been extensively studied in order to understand the var- 
ious anomalous properties of the cuprate superconduc- 
tor. The up-date investigations, however, can not sup- 
ply a definitive evidence to show the stable d-wave su- 
perconducting ground state in these strongly correlated 
modelsi. 

It was known that since the Hubbard model for two- 
dimensions is not exactly solvable, using the Hubbard 
model to study the exchange correlation of the system is 
difficult. On the other hand, the t-J model is weaker 
in studying the long range charge correlation because 
the on-site Coulomb repulsion becomes trivial due to the 
non-doubly occupied projection. A better phenomcno- 
logical model to include stronger correlations is the t-J-U 
modeli. The existence of both J- and U- terms is very 
important in a possible new mechanism of superconduc- 
tivity, gossamer superconductivity, proposed by Laughlin 
recentlji^i^. Both J and U appearing in the model has 
been argued to be the result of the correlations like charge 
transfer processes in the three band Hubbard modeli. 
Dealing with the three-band Hubbard model, however, is 
very complicated and thus, the precise analytical deduc- 
tion from the three band to single band Hubbard models 
lacks. Moreover, an exact reduction from the three band 
model to the single band with both J and U terms are 
more difficult. 

In this paper, we would like to deal with the one-band 
Hubbard model with a hopping asymmetry at finite U. 
In this case, since the on-site Coulomb repulsion is not 



infinite, there may be a fraction of the lattice sites dou- 
bly occupied by electrons. Thus, the U-term is non- 
trivial but can be exactly treated. The difficulty is to 
deal with the kinetic term. We will present a variational 
Monte Carlo calculation for the Hubbard model in a two- 
dimensional square lattice. We examine two types of 
variational wave functions, the resonating valence bond 
(RVB) wave function and anti-ferromagnetic-RVB (AF- 
RVB) coexisted wave function. There are many varia- 
tional wave functions, including AF, RVB and AF-RVB, 
due to different kind of approximation. The early studies 
shows that the results of AF and RVB were contradic- 
tory but the AF-RVB had lower energy- So it is 
quite reasonable to consider the AF-RVB wave function. 
However, the mean field studies prefer RVE^-^^. So we 
included RVB in our work that it might get a compar- 
ison with the mean field studies. It is found that the 
RVB state has a lower variational energy for smaller U 
and tmix which is the hopping amplitude changing the 
local double occupation while the AF-RVB state is fa- 
vorite for larger U and tmix- For both wave functions, we 
check their phase structures. Both wave functions have a 
crossover from metal to insulator states as U /tmix tends 
to infinity. On the other hand, there is an ambiguity to 
find an optimal variational paring parameter A describ- 
ing the RVB feature of the states. In a wide range from 
A = to 1, the variational energies are almost degener- 
ate for the lattice size in our calculation. This leads to 
a difficulty to identify if the metal state is either Fermi 
liquid or d-wave superconducting. 

To understand the physics in the crossover regime more 
clearly, we study an effective model which includes the 
contribution up to the second order of tmix/U. More- 
over, the experience in the t-J model taught us, if there 
is a spin exchange term in the Hamiltonian, the pairing 
variational parameter is much easier to be optimized^. 
For a large U, the spin exchange may be explicitly shown 
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by considering the virtual hopping of electron between 
two single occupied sites. A hopping term changing the 
double-occupancy may be taken as a perturbation as that 
in deducing the t-J model from the single band Hub- 
bard model. In the perturbative deduction from the 
Hubbard model to the t-J model, we see that the J- 
term is actually from the second order perturbation in 
tmix/U by taking the off-diagonal term for the double- 
occupied number in the kinetic term as the perturba- 
tion Hamiltonian. Namely, if the t-term is decomposed 
into Tfjiag + roff-diagj the perturbation process transfers 
Toff-diag to J-term, a virtual hopping process, and only 
Tdiag serves as the real hoppingiS. It can also clearly be 
seen from the canonical transformation deduction of the 
t-J modeili'. In the present case, we can still have a J 
term as in the t-J model while the U term is kept due 
to a non-zero double-occupancy. In the t-J model, due 
to the no-double occupancy constraint, the kinetic term 
is a hopping between a single-occupied site to an empty 
site. For the present model, besides this hopping term, 
a hopping between double- and single-occupied sites and 
a pair hopping between the double-occupied and empty 
sites would be included. Thus, we can derive an effec- 
tive theory at large but finite U, which captures both 
the charge and exchange correlations of the system. We 
can have kinetic, J- and U- terms. However, it should 
not confuse with the t-J-U model mentioned above. In 
t-J-U model, the J is set as a free parameter that means 
J is independent of t and ifi. But in our case, J comes 
from the expansion of tmix/U. In this effective theory, 
there is no a hopping term which changes the local dou- 
ble occupation. It is an extension of the t-J model with 
doubly-occupied sites. 

To work out our model, we use the canonical transfor- 
mation. We find that, to the second order of tmix/U, the 
effective Hamiltonian can be written as the sum over the 
Hamiltonians acting on a subspace of the Hilbert space 
with a fixed double occupied number D. This fixed D 
Hamiltonian including three hopping terms which serve 
the electron hopping from single to empty sites (th term) , 
double to single sites (t^ term), and the paring hopping 
(tp term) ; the U term and J-term as well as various near- 
est neighbor interactions. 

We also study the RVB and AF-RVB variational wave 
functions for this effective model by the variational 
Monte Carlo method. The RVB state is in smaller J 
and U regime while AF-RVB state is favored in larger 
J and U, consistent with the Hubbard model. For both 
wave functions, a first order metal- insulator transition 
may be found^^. Finally, we can plot a phase diagram in 
J-U plane. The regime J = 4:{tmix)'^/U with tmix/U <C 1 
should describe the physics of the crossover regime in the 
Hubbard model. We see that if neglecting the high or- 
der term, this crossover corresponds to a metal-insulator 
phase transition. 

This paper was organized as follows: In Sec. II, the 
detailed deduction of the effective model is provided by 
canonical transformations. In Sec. Ill, the VMC results 



for the original Hubbard model and the effective model 
are presented. In Sec. IV, we give some discussions and 
conclusions. The mean field theory is arranged in the 
appendix for giving some feeling to relate our effective 
model to Laughlin's gossamer superconducting model. 



II. HUBBARD MODEL AND ITS LARGE BUT 
FINITE EFFECTIVE MODEL 

A. One-band Hubbard Model 

We start from the Hubbard model on a two- 
dimensional square lattice where the hopping energy may 
be dependent on the occupation of sites involved^'^. In- 
cluding the on-site Coulomb interaction, this Hubbard 
model reads 

L 

= T + y = T+C/^j.„ (1) 
1=1 

where L is the number of the site; i^i — n.i|7iij^, riia- = 
^l^Cia with Cia a spin-fj electron annihilation operator at 
site i and the kinetic term is given by 

T^Th + Td + T^i^, (2) 
Th = t'ijil - n,s)c\^Cj^{l - Ujs), 

Td = — t^jTligcl^Cjanja-, 

= — t^j nigc\^Cja{l — rija-) 

— tij (1 — n'ig)cl^Cjcrnja-- 

Here T+{T^) creates (destroys) a double-occupied site. 
We assume ijl = = t and t"]^'^ — t™ for the nearest 
neighbor sites and vanish otherwise. 



B. Effective Model 

In large but finite U {U ^ tmix), we can treat the 
Tmix term as perturbation, which leads to the t-J model 
in infinite U limit. An easily pellucid way to arrive at 
the effective model is via a canonical transformation. In 
order to define the canonical transformation, we explain 
our notations. The partial Gutzwiller projection opera- 
tor 

N/2 

Ii{9) = \{{l-{l-gy^) = Y.9°PD=9^, (3) 

i D=Q 
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where < <? < 1 is the Gutzwiller parameters; N is the 
electron number™, D = i^j and 

{ii,...,iD} j 

is a projection operator which projects a state into 
the subspace with a fixed double-occupation number D. 
Pq = n(0) is the full Gutzwiller projection operator and 
n(l) = 1. For convenience, we denote 

Pd{9) = g''PD,Prui9) = E Pd{9). 

D>i 

The first goal of this work is to construct an effective 
Halmitonian H^s and after the partial Gutzwiller projec- 
tion, the projected effective Hamiltonian is given by 

n(5)Heffn(.g) = Y.PDig)H,sPDi9) 

D 

= Y.9'"'PdH,sPd. (4) 

D 

i.e., all the off-diagonal part P]j'{g)HcfiPD{9) — for 
D' ^ D. We shall prove that the result effective Hamil- 
tonian in which all terms keep Z?-invariance is given by 



where 



H,s^Th + Td + Tp + J + V, 



(5) 



where Tp is a pair hopping kinetic energy and J7 is the 
spin exchange as well as various nearest neighbor inter- 
actions, namely. 



T„ = 



J = 



1 

2 



jij(Si • Sj - 



(6) 



For trmx/U < 1, Jy = J « 4t^i^/J7, tp = J. 

The canonical transformation for the Hubbard model 
to the t-J model has been a standard techniquaii. A 
detailed review for the canonical transformation can 
be found in Ref^. Our derivation is a generaliza- 
tion of the D = case. Notice that PoTmixPo' = 
Sd' D±i-PD?mix^Z5±i ^-iid (O remains D invariant, as well 
as U{g)Uig') - H(5.g'), Pno{9)Pnoi9') - Pnoi99')- 
Keeping these in mind, we do a partial projection 
Il{x)HII{x) with X = g^^^ ■ For large N, x is very close 
to 1. A straightforward calculation leads to a rewriting 
of U{x)HIl{x) 



U{x)Hn{x) = Hq{x) + H^pix), 



(7) 



lag (-^ 

D=0 

H^'^Hx) = Pd-i{x)TPd{x)+Pd{x)TPd-i{x).{%) 



The purpose of the canonical transformation is to acquire 
an effecti\ 

PDH[^^Pn 



D^off -^0 = for I? 7^ to the second order of t/U. 

(1) 

off 



an effective Hamiltonian such that PoH^^Pd — 

r(i)p„ - 



This H^^ is defined by 



(1) 



(9) 



As well-knowniiiiS,, S*^^-* is determined by the self- 
consistent condition 

iHl^\x) + [Ho{x),S^^'^]^Q 

and thus the effective Hamiltonian reads 

h'^^ ^ H,{x) + l-[S^'\H^^Hx)] 



(10) 



Solving the self-consistent condition , in a large U 
is given by^StiSi 

PoH^Jf^Po^PoHPo-^PoHP,„HPo, (11) 
P.AxM'JPM « P,n{x')HP^,{x') 

+ ^P,,ix')HP,HP,,,{x').{12) 

The approximation 'w' in pi|) and (|12|l means the exact- 
ness is up to the second order of t/U. Namely, the third 
term in pU|) has been neglected. In fact, the off-diagonal 

part PoHj^g Pri^ vanishes also only up to the second order: 
r(i) 



PoKrP^, = PoHP„,HPoHP,, 
1 
U 



^T_r+r_ - o{tj) (13) 



is of the third order. The second terms of pi|l and p2f) 
may be calculated and given by 



U 



PaHP,j,HPo 



1 

U 



PqT-T^Pq ~ PaJPoi 



D=2 



— Prj^HPoHPjj^ - —Pjj^T+T^P,^^ w P^jTpP^^ (14) 

Thus, up to the second order, we have 

PoH^gPo - PoH^sPo, (15) 
Pr,, (x)i7«P,, {x) « P,, {x'){H + Tp)P,, (x2).(16) 

where the approximation 'w', besides up to the second 
order, also means the three and more sites processes are 
neglected. 
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If the non-double occupied constraint is imposed, H16() 
vanishes because it is related to the double occupation. 
Ea. p5|l gives rise to the common t-J model. How- 
ever, if the double occupation is allowed, we have to 
deal with (|16() . In fact, one can repeats the canoni- 
cal transformation to H16(l . We would like to require 

(2) 

an effective Hamiltonian H^g whose off-diagonal part 

PiH^'^^Pd = PdH'§Pi = 0{tJ) for £) > 1. For this 
purpose, one writes 

Ii{x)H'^J^Il{x) = PoH'^^Po + Ho{x^) + (17) 
where 

Mx^) = Pi{x^){H+Tp)Pi{x^)+P^,{x^){H+Tp)P^,{x^). 
We do a canonical transformation and define 

i^i^) = e^-'^'^(.)i^«^(x)e-^'^ (i8) 

where 5^^^ is required to satisfy Po<5''^' = •S'-^^Pq — such 
that PoIi{x)H^^I{{x)Pii is invariant under the transfor- 
mation and it is sclf-consistently determined by 

iH(^\x') + [H,{x^),S(^'>]^Q. 

Hence, similar to (|10ll . one has 



..(19) 



Projecting H^'^ — H(x)i?^g H(x) and repeating the sim- 
ilar procedure to deduce (|15|l and (|15|l . one arrives at 

PoH^l^Po ~ PqH.sPo. 



r(2) 



Pi{x)H^l^Pi{x) 

r(2) 



x'PiHasPi, 



PvMKl^Pmi^) « PvA^niH + 7;)P„,(x^), 

for a large U, where the three site processes have been 
ignored. 

Repeating this procedure, we finally have 



n{x)Hli ^n(x) 



D=0 



PdH.sPd 



H(g)i7effH(g). 



(20) 



The last equality is because i/eff is i^-invariant. The 
Gutzwiller parameter is g but not x because we are doing 
the partial projection in each time canonical transforma- 
tion. Thus, we end the proof of I0J and ||SJ). Moreover, 
we see that, in a partial Gutzwiller projection, the vari- 
ational ground state energy is given by a polynomial of 
the Guztwiller parameter g in power of 2D. The coeffi- 
cient of g^^-term is the ground state energy of the system 
with a fixed D. Using g as a variational parameter may 
be convenient for the numerical simulations. In the orig- 
inal Hubbard model, the change of the double occupied 
number is allowed. We see here that the allowance of 
this change in a large U is very small. After neglect- 
ing the three and more sites processes, the probability of 
the change of D is in the third order of t/U as eq. (|13f) 
shown. Considering the fixed D processes may be helpful 
to numerical simulations. 



III. VARIATIONAL MONTE CARLO RESULTS 
A. Variational Wave Functions 

The variational wave functions we would like to study 
are so-called the partially projected RVB state = 
Pd\BCS) and the partially projected AF-RVB state 
= Pd\AF-BCS). The BCS state is defined by 

\BCS)^\{(uu + vucyi^^m, (21) 

k 

where Uk and follows the standard BCS form 

a{k) = ^ = -4V' (22) 

Uk ik + Ek 



= -2(cos + cos ky) - fi, Ek = \ 



for the d-wave pairing parameter A^. = A (cos fc^;— cos ky). 
The AF-BCS coexisted state \AF-BCS) is defined hyi^ 



\AF-BCS{A,,Aaf,fi))=l[{u^^^ 



k.s 



(X exp 
where 



|0), 



(23) 



4"^ 



,(±) 



±Ad7fc 



and Ek = ^J^l+^f, = -2{coakx + cos ky) and 7fc 
2(cosfc2; — COS ky); and 



with 



4ct''^ ^ O^kaCAka- - l^kcrCBka, 
d[ = PkaCAka + ^kaCBka, 



OLka 



aka 



(25) 




(26) 



CAk(j{cBk(j) is the electron operator on sublattice A(B). 



B. Hubbard Model with Asymmetric Hopping 

We first make a variational calculation for the original 
Hubbard model. The energy we want to minimize is given 
by 

^ C/d+^ii^^^^^^^#±^ (27) 
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where y = g^, Njj = {iPd\iJ-'d) for the partially pro- 
jected RVB state \tpD) = Pd\BCS) or the partially pro- 
jected AF-RVB state |Vd> = Pd\AF - BCS). 

The average double occupation number d is given by 



And 



T, 



h(d),D 



rjiraix 
^D,D + 1 



{'4'D\Th{d)\'fpD) 

Nd 

{'(pD + l\T+\i}D) 



(28) 



D+l) 



D+1 



(29) 



D+l 



Let {!«£))} be a set of the basis in the configuration 
space with a fixed D. The normal factors Nd is given by 



an 



(30) 



where Aaj^ is just the determinant of the configuration 
aD- We are not able to calculate Njj exactly. We use 
the approximation by taking all probabilities to 
be the same22iSi. Thus, at half-filling 



[{N/2 ^ Dy.]^Dl{L - N + D)\ 



(31) 



for the RVB case. In the AF-RVB case, the lattice is 
divided into two sublattices A and B respectively and 



N 



AF-RVB 



- E 

Na-jNaiNadNae 

{L/2y. {L/2y. 

Na^ \Nai\Nad^.Nae^- Ng-t \Nbi INbd^Nbe ! 



(32) 



where the configurations {NA■\^NAl^NAD^NAE',NB■\^ 
Nbi,Nbd,Nbe) are corresponding to numbers of spin- 
up, spin-down, double occupancy and empty sites for 
each sublattice and subjected to the following con- 
straints: 



Na]- H 


^Nai^ 


-Nad-^ 


'Nae - 


N 




-Nbi ~ 


~Nbd^ 


'Nbe = 


N 




+ iVsT 


^Nai- 


^Nbi - 


N 


Nad^ 


^Nbd 


^Nae- 


^Nbe - 


D 



(33) 



By using the variational Monte Carlo methodic, we 
calculate the variational energy H27(l by optimizing the 
variational parameter A. The term Ud is not dependent 
on A. For the projected RVB wave function, Th.o-Td.o 
and T™^j^^ for several D are depicted in Fig. ^ (The 



energy unit t = 1 is used in all figures through the pa- 
per.) The lattice sizes are 10 x 10 and 12 x 12, respec- 
tively. We use periodic-antiperiodic boundary condition 
to avoid the degeneracy in Brillouin zone. All data are 
calculated with more than 10"* Monte Carlo samples. Al- 
though there is a minimum in T^jj^^j^i around logj^g^ — 0: 
the total kinetic energy is minimized after logj^pA < —1.0 
because the minima oiT^j^^ jy are in after logj^g^ < —1.0. 
Unfortunately, we see that there is very broad minimal 
flat in variational energy from A = to logj^gA ^ —1.0. 
Thus, we can not distinguish the metal state is either the 
Fermi liquid or superconducting state. For the AF-RVB 
wave function, the trend of Th^n + T^^d and T„iix,D,D+i 
is different. However, in the total energy 123 T^.d + Td.o 
dominate. So the situation is like the RVB case. The 
parameters we use are logj^Q A = logj^Q A^ = —0.6 and 
the optimal log^o Aq/ w — O.e^iSS. 



Compare the variational energies of the two wave func- 
tions, we find that the system is in the projected RVB 
state for small U and tmix while it is in the projected AF- 
RVB coexisted state for larger ones. Tablc^ shows the 
transition when Tmix = 0.6. The critical line is shown in 
Fig. El 

To understand the phase diagram of the system, we 
shall calculate the optimal average double occupied num- 
ber d for an appropriate wave function (RVB or AF-RVB) 
for given tmix and U in the optimal parameters (A or A^ 
and Aaf). Substituting into l|?7jl and ehminating y, 
we get the function E(d). Then, identifying the minimum 
of E over d, we get the optimal do and Eq. If do = 0, the 
system is in insulating state while if do > 0, the system 
is in metal state. There is a second order phase metal- 
insulator transition in tmix = as show by Fig. ^t,. The 
critical interaction Uc{0) is spotted in Fig. |2] However, 
when tmix > 0, 



dEjd ) , 
dd 



= u- 



(ri-ro) + i-^T™^ 



(34) 



where C is a constant. For any finite tmix, no matter 
how large U is, it can be found there exists a do > 
so that = 0. As instances, in Fig. |SIb)(c), we 

plot the d-E curves for tmix — 0.8 and U — 10 for the 
RVB state (Fig. E^b)) and the AF-RVB state (Fig. O^c)). 
The dashed curve in Fig. El gives the values of {tmix , U) 
where the do = 0.01. For a sufficient small do, the system 
becomes a practical insulator and therefore, there is a 
metal-insulator crossover as showed by the shade area in 
Fig. El Due to small U, the RVB region is in metal phase. 
The AF-RVB region is divided into two phases. For a 
given J7, the system is in the insulating phase when t^ix 
is small enough while in metal state when tmix is large. 
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u 





1 


2 


3 


4 


5 


6 


4.316 


Ervb 
Eaf-rvb 


-0.835 
-0.795 


-0.603 
-0.564 


-0.405 
-0.373 


-0.242 
-0.221 


-0.115 
-0.109 


-0.0245 
-0.0425 


-0.00202 
-0.0205 


-0.8324 
-0.8324 



TABLE I; The transition of Hubbard model between RVB and AF-RVB at tmi^=0.6. The critical U is 4.316. 



C. Effective model 

We now begin to examine the effective model. The 
energy we want to minimize is defined by 

where d defined by H28() . 

For fixed D,Nd defined by (EOJ), (EU and And 

Th,D, Td^D, Tp^D, Jd is defined by 



T, 



h{d,p),D 



{iw\Th(d.,p)\^D) 



N 



D 



Jd 



{^d\J\^l 
Nd 



(36) 



Our strategy is that using the variational Monte Carlo 
method to minimize Ed for fixed D and fixed electron 
number N at the half-filling by varying the variational 
parameter log^pA. Then, draw the curves E as the func- 
tion of d through eqs. and (^5)) . to read out crit- 
ical Uc and dc from the shape of the curve for differ- 
ent model parameters U/t and J/t. At the moment, al- 
though we still use J = At'^^^/U, we do not not restrict 
at tmix/U <C 1. The comparison to the Hubbard model 
is only valid in the region tmix/U <^ 1. 

Our variational Monte Carlo carries out on square lat- 
tices as in the Hubbard model above, with sites L from 
10 X 10 to 12 X 12. A periodic-antiperiodic boundary con- 
dition is used. All data are calculated with more than 
lO"' Monte Carlo samples. In the half-filling, we set the 
chemical potential — 0. The ground state energies 
Ejy are calculated. We show Si ■ Sj for D=0 and D=l 
varying as log^gA in Fig^t. for the RVB state. The no- 
double occupant energy D=0 is the variational ground 
state energy of the common t-J model. Our result is well 
consistent with the known resultsSiii. We calculate Ed 
up to the largest D = L/2 — 1, and find that all these 
energies are almost degeneracy in wide range between 
—0.5 < logj^pA < 0.0 . Using the Monte Carlo estimat- 
ing energy Ed on 10 x 10 lattice, we approximate E in 
by finite sum for £> = 49 and logipA = -0.5. The 
error bars for independent Monte Carlo initial configura- 
tion are in order of 1% and we do not show them. 

The energy of the AF-RVB wave function also can be 
calculated by variational Monte Carlo method with op- 
timizing both of the model parameters A^/ and A^.We 
show Jd for D— and D=5 in Fig^ and Fig^. The re- 
sults of D=0 corresponding t-J model at half filling. Our 
results are consistent with the known results. One can 



see that for D=0, the energy minimum locates in a deep 
valley. 

We also analyze the two wave functions' finite-size scal- 
ing of D=0 which corresponds to Heisenberg model. The 
results are show in Fig|Sl All the data but the 16x16 
of AF-RVB, which is only one datum since it is required 
very long time to get one result, are average of 5 inde- 
pendent calculations. One can see that for Heisenberg 
model the energies of AF-RVB are deeper than those of 
RVB. 

For a pair of fixed J and U, we can compare the 
variational energies corresponding to both wave func- 
tionsfTab HH) . In this way, the J-U plane can be divided 
into two regions: RVB and AF-RVB, which is similar to 
the case in the last subsection for the Hubbard model. 

For a given type wave function, we look for the pos- 
sible metal-insulator transition. First, like the case in 
Hubbard model, there is a second order phase transi- 
tion when J = 0. If J > 0, due to the vanishing of 
rprmx iQYui, there are first order phase transitions in a 
given type wave function. Fig. (Ogives an example of the 
first order metal-insulator transition. In Fig. |S1 we show 
the relation between the critical Uc or Jc and the critical 
double occupied concentration dc for the RVB wave func- 
tion. In this way, we can determine the critical Jc — Uc 
line in J-U plane. Figs. I3a)(b) show the critical Jc — Uc 
lines for the RVB and AF-RVB wave functions. They are 
quite similar. Combining these two phase diagrams to- 
gether with the region-dividing picture mentioned above, 
we depict the comprehensive phase diagram (Fig. Eta)). 
In the RVB region, due to small J and U, the system 
is in metal state. In the AF-RVB region, the system is 
basically in an insulating phase. For J > 0.5, there is a 
phase which may be a AF-RVB metal state. Since the 
optimal variational parameters A and A^; are not zero, 
the metal state may be a superconducting state. Con- 
verting J t,nix (see Fig. Elb)), we find that for small 
trnix, the phase diagram is consistent with the crossover 
picture in the Hubbard model . 



IV. 



DISCUSSIONS AND CONCLUSIONS 



We have investigated the Hubbard model with the hop- 
ping asymmetry and deduced an effective theory for large 
but finite U. Based on two types of the variational wave 
functions, the phase diagram of both models are depicted 
by the variational Monte Carlo method. For the Hub- 
bard model, we found it is difficult to determine the ex- 
act critical boundary of the phase transition of metal- 
insulator. Moreover, the superconducting behavior in 



7 



u 


1 


2 


3 


4 


5 


6 


3.512 


Ervb 


-0.656 


-0.491 


-0.381 


-0.34 


-0.34 


-0.34 


-0.3494 


Eaf-rvb 


-0.616 


-0.457 


-0.365 


-0.34943 


-0.34943 


-0.34943 


-0.3494 



TABLE 11: The transition of t-J-U model between RVB and AF-RVB at J=0.3. The critical U is 3.512. The critical U of RVB 
M-I transition is 3.8, and that of AF-RVB is 3.45. 



the metal phase was not clear. The effective model is 
a finite but large U extension of the t-J model. This 
model captures both the charge and exchange correla- 
tion. The phase diagram of this model clearly shows a 
metal-insulator phase transition. Due to non-zero opti- 
mal A and A^, the metal state may be superconducting, 
which leads to the possibility of the gossamer supercon- 
ductivity in the framework of the hopping asymmetry 
Hubbard model. 

The relation to the gossamer superconductivity can 
also be seen from the mean field state of our theory. The 
basic idea to go this mean field state has been explained 
in our previous preprin'feS^. Here we present a renewed 
version of the mean field state. We only try to show our 
mean field theory may formally be equivalent to Laugh- 
lin's gossamer superconducting model. We do not intent 
to go more analysis such as the stability of our mean field 
state against other possible instabilities before we work 
out some more sophisticated issues. We put this formal 
identification into Appendix A. 
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X (-l)"(Aj,c..Cj5 + A,jc]^cJJ 

+ U^Ui^nn -^2Jij{l ~ A)ni^nii (Al) 

- X J,,(A|A,,f + i(l - B)\x^,\')in, + n,) 



where the parameters are given by 



7(1) 



,(2) 



(2) 



B 



- Jij , 



t%+ti^-{l-B)x,^J^, 



(A2) 



A and B are the variantional parameters to be deter- 
mined. On the other hand, we write down Laughlin's 
gossamer superconducting Hamiltonian 



(A3) 



APPENDIX A: MEAN FIELD STATE 



We outline the mean field state of our model in this 
appendix. Due to the paring hopping is of the order J, 
we neglected it in our mean field theory. Introducing 
two correlation functions Ay = (ci|Cjf — Ci|Cj|)o, Xij = 
('-It'^jT +'^Ix'^ji)oi the U{1) symmetry of Hcs is broken by 
a decomposition of the four particle terms^'^. According 
to Ay and Xij : the mean field Hamiltonian of |SJ| is given 



where is renormalized chemical potential, Ek = 
Vi^k - Mfl)^ + A^ and bi^„ ^ U{g)bk„U-^{g) for 6kT 
MkCkT + i'kclj.| and 6ki — WkCkj, — ''^'k'^^-k^ annihilate the 
BCS state. Expressing exphcitly (|A3|I by the electron 
operatorsSiS^, we have 

- a(3nignjg]{-iy {AljCaCjs + Ayc]^4^) 

-f /JG^n^in,! -mgA^ (A4) 
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where a = I — g and /3 = {1 — g)/g and 



fe 

k 

= ^^Ukz;ke'''-^"^'-'"^\ (A5) 

k 

and Ug = i^Ei.Ek[{2P + P')ul + {2a^a')vl],^^G = 
■W Ek^fe[(2a + l>fc-"fc]- If we identify the t-J-U model 
to the gossamer superconducting model in the mean field 
level, one requires 

A = aP, B = 2a/3, 

tGW = -t%-{\-aP)x,^J^„ 

= tij + 4 - ( 1 - 2«/3)x,^ J^J , ( A6) 



and HR + ^lG = J(12A|A,|2 + 8(1 - B)|xrP) + ^^..U = 
Uc + SJil-A). 



Although we have made a formal equivalence between 
our mean field state Hamiltonian to Laughlin gossamer 
superconducting Hamiltonian, we note that the hopping 
parameters t'^^^^^) have run out of the practical range in 
the real materials. Thus, to show the system described by 
the t-J-U model has a gossamer superconducting phase 
described by Laughlin gossamer superconducting Hami- 
tonian, a renormalization group analysis is required. We 
do not touch this aspect in this work. However, we can 
believe there is such a superconducting phase in our the- 
ory if J7 < Uc because the superconducting paring pa- 
rameter is determined by the optimal exchange energy 
as in the common t-J model. The renormalization of 
the hopping parameters is believed to affect the normal 
dissipation process only. 
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/o, 



(d) 



FIG. 1: (a)(b)The kinetic energies Th,D + Td,_D and 
D.D+i of RVB as functions of logj^gA. (a) r-,5,6 +T'+,6,5 
(squares), rh,5 + Td.s (circles) and Th,6 + Td.e (triangles) in 
a 10 X 10 lattice, (b) T_,7,8 + 7+, 8.7 (triangles), Th,7 + Id,? 
(squares) and Th,s + Td,s (circles) in a 12x 12 lattice. (c)(d)Tlie 
kinetic energies Th^n + Td,D and Tmix,D,D+i of AF-RVB as 
functions oflogjoAd andlogj^pAa/ inalOxlO lattice. (c)Th_5 + 
7d,5.(d) r_,5,6 + r+,6,5. Note that for the different trend be- 
tween Th,D + Td,D and Tmix,D,D+i, we change the view of 

Tmix,D.D + l- 
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FIG. 2: The possible phase diagram of the Hubbard model. 
The solid line divides the tmix-U plane into the RVB and 
AF-RVB regions. The RVB region is in metal state (RVB- 
M). There is a crossover from metal to insulator in AF-RVB 
region( AF-RVB-M to AF-RVB-I). The spot on the U-axis 
{tmix = 0) is the metal-insulator phase transition point (7c (0). 
The shade area is the crossover region and along the dashed 
curve, do ~ 0.01. 
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FIG. 3: The d-E curves of Hubbard model, (a) tmix = 0. 
The different curves are corresponding to different U. Lower 
curve has a smaller U. The second phase transition happens 
at Uc{0) = 4.5678. (b) and (c) are drE curves for RVB state 
and AF-RVB state at tmix = 0.8 and U = 10, respectively. 
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(c) 

FIG. 4: The variational energy for fixed D varying as the vari- 
ational parameter logj^pA for RVB(a), logioAd and logioAaf 
for AF-RVB(b)(c). 




FIG. 5: Finite-size scaling of RVB at logio A=-0.5(a) and AF- 
RVB at logioAd=-0.6 and logioAa/=-0.6(b).Where L is the 
lattices size. 




FIG. 6: The total energy E averaged in RVB variational wave 

function varying with the double occupant concentrate d for 
J = 0.3 and (/=3.804, 3.8007, 3.7974 from the upmost to the 
lowest, respectively. The critical Uc = 3.8007 and dc = 0.0135 
in this case. 
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FIG. 7: The critical concentration dc for the RVB state (a) 
Uc-dc curve; (b) Jc-dc curve. 
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FIG. 8: Phase diagrams for different variational wave func- 
tions, (a) for the RVB wave function; (b) for the AF-RVB 
wave function 
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FIG. 9: Phase diagrams of the effective model (a) the phase 
diagrams of the system in the J-U plane. (b) converting (a) to 



